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The issue of intermittency in numerical solutions of the 3D Navier-Stokes equations on a periodic 
box [0, L] 3 is addressed through four sets of numerical simulations that calculate a new set of 
variables defined by D m (t) — (vjQ 1 Q r n)° Lm for 1 < m < oo where a m — 4m m L 3 and [Q m (t)] 2m = 
L~ 3 f v \u;\ 2rn dV with vjq = vL~ 2 . All four simulations unexpectedly show that the D m are ordered 
for m = 1 , 9 such that D m +\ < Dm- Moreover, the Dm squeeze together such that D m +\/ D m /* 
1 as m increases. The first simulation is of very anisotropic decaying turbulence ; the second and third 
are of decaying isotropic turbulence from random initial conditions and forced isotropic turbulence at 
constant Grashof number respectively ; the fourth is of very high Reynolds number forced, stationary 
isotropic turbulence at up to resolutions of 4096 3 . 

I. INTRODUCTION 
A. Background 

Intermittency in both the vorticity and strain fields is a dominant feature of developing and devel- 
oped turbulence. It has been studied extensively both experimentally (Sreenivasan 1985, Meneveau and 
Sreenivasan 1991) and numerically (Kerr 1985, Jimenez, Wray, Saffman and Rogallo 1993, Donzis, Yeung 
and Sreenivasan 2008, Ishihara, Gotoh and Kaneda 2009, Donzis and Yeung 2010, Donzis, Sreenivasan 
and Yeung 2012, Yeung, Donzis and Sreenivasan 2012). Statistical physicists generally use velocity struc- 
ture functions to study this phenomenon and have diagnosed the degree of intermittency by how much 
the velocity structure function exponents ( p differ from linear when p > 3 (Frisch 1995, Schumacher, 
Yakhot and Sreenivasan 2007, Boffetta, Mazzino and Vulpiani 2008, Pandit, Perlekar and Ray 2009). 
The standard way to quantify equal-time, multi-scaling exponents is a challenging numerical task (Ar- 
neodo et al. 2008, Ray, Mitra and Pandit 2008, Ray, Mitra, Perlekar and Pandit 2011). The multi-scaling 
approach is even more challenging for the 3D Navier-Stokes equations 

d t u + u\7u = i/Aw-VP divu = (1) 

because the velocity field u(x, t) and pressure P(cc, t) evolve in time, so, in general, time-dependent 
structure functions must be used to study dynamic multi-scaling (Ray, Mitra and Pandit 2008, Ray, Mitra, 
Perlekar and Pandit 2011). This paper will introduce an analysis of some new and existing numerical 
computations that gives new insights into and provides a new method for distinguishing alternative 
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regimes of behaviour in the vorticity field. To explain the nature of these regimes, let us consider the 
vorticity field uj = curlu on a finite periodic domain [0, L] 3 within the setting of the volume integrals 
which define a set of frequencies 

/ r \ 1 / 2m 

n m (t) = I L~ 3 J \u\ 2m dv) , 1 < m < oo . (2) 

Some recent work has centred around a dimensionless re-scaling of the Q m such that (Gibbon 
2010, 2011, 2012a,b) 

D rn (t) = (zuQ 1 Q m ) am , a rn = - — ni — , (3) 

4m — o 

where is a fixed frequency defined by = vL~ 2 . The origin of this re-scaling, valid for both the 
Navier-Stokes and Euler equations, has been explained elsewhere (Gibbon 2011, 2012a,b) where it has 
been shown that, with additive L 2 -forcing, weak solutions obey the time average up to time T 

{D m ) Umeav .<cRe^O(T-') . (4) 

The first in the hierarchy, D\ = Wq 2 Z, is proportional to the global enstrophy Z = Q 2 and may be 
insensitive to deep fine-scale fluctuations. The higher D m may be more sensitive so their measurement 
over a wide range of m could be a useful diagnostic of intermittency. However, the end of the sequence, 
Doo^), is hard to measure numerically, especially in highly intermittent flows. While Holder's inequality 
enforces a natural ordering on the frequencies ^ m such that Q m < Vt m +i for 1 < m < oo, no such natural 
ordering is enforced on the D m because the a m decrease with m. Thus there are two possible regimes : 

D m +i(t) < D m (t) , (regime I) , D m (t) < D m +i(t) , (regime II) . (5) 

The issues to be addressed in this paper in our four numerical simulations of the 3D Navier-Stokes 
equations are : 

1. Which of these regimes is favoured or is there an oscillation between them? If one regime is favoured, 
are the D m well separated? What is the role of the enstrophy D\l 

2. Are these processes m-dependent? 

3. Are they i?e-dependent? 

4. Are they dependent upon initial conditions? 



B. Simulations used for tests 



An important point with respect to numerical simulations of the scaled higher order moments D m (t) 
is that their ratios might converge better than their actual values. This is consistent with the results 
reported in Yeung et al. (2012) and Donzis et al. (2012) where convergence for the ratios of higher-order 
vorticity and dissipation (strain) moments were obtained, even when the statistics of the individual 
moments showed no evidence of convergence (Kerr 2012a). This answered a problem first raised in Kerr 
(1985) where it was noted that in forced simulations at modestly high Reynolds numbers, the averages of 
the vorticity and strain moments above 6-th order did not converge. The determination of the D m (t) in 
simulations is not difficult whereas, in contrast, traditional numerical tools such as higher-order structure 
functions require a combination of larger domains and finer resolution than is currently feasible. This 
paper will calculate and compare the D m (t) from four data sets: two where the average kinetic energy 
E = L~ 3 j v ^\u\ 2 dV decays in time, and two where E is held approximately constant by forcing at the 
low wavenumbers. The first is a unique data set from a computation in which fully-developed turbulence 
forms from the interaction of two anti-parallel vortices and whose kinetic energy E decays strongly 
after the first peak in the normalised enstrophy production — S u . Because this calculation has not been 
fully described before, some introductory discussion is provided at the start of that section (^H|. The 
other three data sets represent more traditional decaying and forced homogeneous, isotropic numerical 
turbulence. In the decaying calculations in 3ll| and the decaying and forced calculations in £111 the 
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moments have been determined relatively continuo usly in time which makes a helpful comparison with 
the results of ^ For the fourth data (4096 3 ) set of pv| (Yeung et al. 2012, Donzis et al. 2012), a similar 
conclusion is reached by studying the dependence of the average value of D m on the Reynolds number. 

An advantage of the first data set described in jjll]is that the predicted convergence properties of ratios 
of the D m (t) can be tested for a calculation with huge fluctuations in the production of enstrophy, and 
therefore in the higher D m (t). That the calculation eventually exhibits traditional turbulent statistics 
and spectra is a bonus in justifying its use. However, this new initial condition is very specialized and 
any trends need to be confirmed using a more traditional decaying homogeneous, isotropic data set, 



which is the purpose of the second calculation discussed in { III §III| also contains forced simulations of 



homogeneous and isotropic turbulence at constant Grashof number. Finally, the fourth calculation in pVJ 
provides validation from a forced, massively parallel, pseudo-spectral calculation (4096 3 with R\ w 1000) 
calculation to show that these trends are not restricted to low or moderate Reynolds numbers. Assessing 
the scaling of moments of intermittent quantities such as vorticity, strain rates or velocity gradients has 
been a critical component of characterizing and understanding intermittency. Of particular interest is 
how these moments scale with the Reynolds number, which is typically high in applications. At the same 
time, different orders provide information about fluctuations of different intensities. Low and high-order 
moments, for example, are associated with weak and strong fluctuations. Thus, the understanding of the 
scaling of the moments D m , especially at high m, can also shed light on the nature of intermittency and 
the most extreme events in turbulence. 



C. A summary of results 



The simulations described in £jll| £111 and §IV[ and illustrated in Figs. [T][2]and|5j each observe that a 



strict ordering of the D m occurs, as in regime I ; namely D m +i < D m (on log-linear plots). To assess the 
significance of this, we write down the relation D m +i < D m in terms of Q m and use Holder's inequality 
< ^m+i on the extreme left hand side 

vj^Qm < w^Qm+i < (tz7 " 1 ^ m ) arn/ar " +1 . (6) 

As m — » oo, a m \ a m +i, and so (|6| shows that Q m +i/Q m \ 1. Thus, in regime I the Q m must be 
squeezed together for high m. In terms of the D m (|6| is written as 

D^/ a ^<D m+1 <D m . (7) 
While respecting the ordering D m +i < D m , D m +i is squeezed up close to D m as m — >• oo 

lim %±Vl. (8) 

This squeezing phemonenon is observed in all four data sets where the D m -curves lie very close for m > 3 
as in Figs, [I] [2] and [5] Moreover, the values of D\ in all four simulations lie far above the rest of 
the Dm giving rise to a suggestion, explored in that a depletion of nonlinearity is occurring which 
could be the cause of Navier- Stokes regularity. The most extreme intermittent events are represented by 
moments at increasingly large m. Our results suggest the saturation of these high order moments. This 
is significant as it constrains the shape of the tails of the PDF of vorticity which has been the focus of 
intense investigations (Kerr 1985, Jimenez et al. 1993, Donzis et al. 2008, Ishihara et al. 2009, Donzis 
and Yeung 2010, Yeung et al. 2012, Donzis et al. 2012 ). The fourth data set (forced, stationary, isotropic 
turbulence), the results of which are displayed in j |IV| furnishes us with the opportunity to compare these 
results with other results on intermittency available in the literature. For example, within the multifractal 
model, Nelkin (1990) found that normalized moments of velocity gradients scale as 

«) / (ulf 2 ~ Re x d » , (9) 

where d p is obtained from the multifractal spectrum and (•) is the usual notation for the statistical average. 
Using the well-known result (v%.) ~ (Uo/L) 2 Re\ due to the dissipative anomaly, it is readily shown that 
(v%.) ~ Re\ pJrdp . Our interest lies in the limit p — » oo where it can be shown that lim^oo d p /p = c. 
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FIG. 1: First panel: The initial condition is characterized by long anti-parallel vortices with a localized perturbation for the 
Re = 4000 reconnection calculation. Second panel : Plots of the kinetic energy decay E, the enstrophy Z and its production, 
normalised to be consistent with experimental measurements of the velocity derivative skewness — S u . Z grows until t « 270, 
while E is always decaying. Third panel : Curves range from m — 2 to m — 9 and include the maximum of vorticity . The 
normalized enstrophy, D±, is far above the log-scale given here, so it is omitted. The Dm are ordered for all values of Re and all 
times. 



The constant c is given by c = 3(1 — X><x>)/(3 + £>oo) with representing the limit lim^ 
generalized dimensions V q (Nelkin 1990, Hentschel and Procaccia 1983). Clearly moments of the form 



T> q of the 



(uP) 1 ^ saturate at high p, consistent with (|8j). While experimentally it is difficult to measure reliably, 
its value appears to be smaller than 1.0 (Meneveau and Sreenivasan 1991). The ratio of successive orders 
is 



« +1 f (p+1) / «y /p 



Re 



(i+d p /p)-(i+d p+1 /( P +i)) 



(10) 



The limiting behavior of d p shows that lim p ^ oc [(l + d p /p) — (1 + d p +i/(p + 1))] = 0, and therefore the 



ratio on the left hand side of (10) tends to a constant independent of p and Re\. This is consistent with 
the squeezing together of the Vt m and D m . 



II. THE FIRST SET OF SIMULATIONS 



The new vortex reconnection calculation displayed in this section addresses the following long-standing 
numerical question : Can an initial condition with only a few vortices generate and sustain fully-developed 
turbulence through reconnection events in a manner similar to how turbulence forms in aircraft wakes or 
when anti-parallel quantum vortex lines reconnect numerically (Kerr 2011)? Because this initial condition 
has not been fully explained before, some of its unique features are now described. The three directions 
in the flow are : (i) each initial vortex primarily points in the indirection ; (ii) separation between the 
vortices lies in the z-direction ; (iii) they propagate in the x-direction. Due to the anisotropy of the 
flow, an anisotropic mesh and domain are used with a L x x L y x L z = 2ir(2 x 8 x 1) domain and a 
n x x n y x n z = 512 x 2048 x 512 mesh, plus symmetries used in the y and z directions. The quantum 
vortex work in Kerr (2011) has shown that the two most important properties are that their initial 
perturbations need to be localized far from the periodic boundaries and their initial profile and direction 
should be balanced so that they are neither internally unstable nor prone to the shedding of waves or 
vortex sheets. The method grew from addressing calculations identified by Bustamante and Kerr (2008), 
where a single sign of the initial local vorticity had not been imposed rigorously. To accomplish these 
goals, the following four changes have been made to the initial condition in Bustamante and Kerr (2008) : 
i) The perturbation is strongly localised near the symmetry plane using the trajectory given in Kerr 
(2011) ; ii) The initial vorticity profile uses the solution of a two-dimensional vortex with a smoothed 
core ; iii) The direction of vorticity chosen at grid points follows the path of the nearest point in 3D space 
on the prescribed trajectory of the central vortex line ; iv) The vortices need to be twice as long as in 
any previous anti-parallel study. Fig. [I] (left) shows the vorticity as the initial instability saturates. By 
using these choices, regions of negative vorticity and vortex sheets on the y = plane, as described by 
Bustamante and Kerr (2008), are eliminated. 

The ultimate goal is to generate turbulence with a persistent -5/3 energy spectrum and additional 

turbulent statistics, including the experimental velocity derivative skewness S u = (u x ) / {u x ^f > ^ 2 \ which 
is equivalent to the numerical normalized enstrophy production. The latest infinite .Re-estimates of S u 
from forced turbulence calculations (Ishihara et al. 2009) find — S u w 0.68, consistent with experimental 
values of — S u ~ 0.5 — 0.7. Early numerical calculations showed that the S u tended to overshoot the 
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FIG. 2: Plots versus time t of the total kinetic energy (first panel, black curve), the enstrophy Z (first panel, blue curve), the 
normalized enstrophy-production rate — S u (first panel, red curve), D m for 2 < m < 9 (second panel, blue to brown curves), and 
Doo (second panel, dark green curve) for our DNS of decaying, 3D Navier-Stokes turbulence ; the value of D\ is very high, so it is 
omitted. The third panel is of statistically steady forced turbulence at constant Grashof number. The mean values of D m in the 
statistically steady state are as follows: (£>i) = 3.1 X 10 11 , (D 2 ) = 5.5 X 10 4 , (D 3 ) = 1.1 x 10 4 , (D 4 ) = 6.6 X 10 3 , (D 5 ) = 5.2 x 10 3 , 
(D 6 ) = 4.5 X 10 3 , (D 7 ) = 4.1 x 10 3 , (D 8 ) = 3.9 X 10 3 , (D 9 ) = 3.7 X 10 3 , and (D^) = 3.0 X 10 3 . Zooming in to the right panel 
makes it clear that D m |i < D m for all values of m considered. 



early experimental values of — S u ~ 0.4 — 0.5 before settling to the expected value (Orszag and Patterson 
1972). The second panel of Fig. [I] confirms this trend for the anti-parallel calculation with — S u first 
rising to — S u — 1.9 at t ~ 45, then falling abruptly to — S u ~ 0.6 at t ~ 80, continuing to fluctuate 
strongly between 0.6 and 1.2 for 100 < t < 250, and finally decaying to —S u « 0.6. The full details, 
plus the relationship between the variations in D m in Fig. [I] and the development of swirling, turbulent 
vortex rings, is the topic of another paper (Kerr 2012b). 

In Fig. [I] note that for all times, all of the lower order D m (m = 1, . . . , 9) bound each higher-order 
D m (on a log-scale) which can be expressed as D m +i(t) < D m (t), thus favouring regime I as in The 
enstrophy Di lies far above all of the other curves and has been omitted. Next note a strong increase in 
the growth of each of the D m , including D^, up until t « 16. This is the period when this calculation 
has nearly Euler dynamics, where the effects of viscosity compared to nonlinear growth are minimal. The 
growth of the D m (t) in true Euler dynamics is the topic of another paper (Kerr 2012c). 



III. THE SECOND AND THIRD SET OF SIMULATIONS : DNS RESULTS FOR 
HOMOGENEOUS, ISOTROPIC TURBULENCE 

Data from two direct numerical simulations (DNSs) of homogeneous, isotropic 3D Navier-Stokes turbu- 
lence is now presented. Both of these simulations use a pseudospectral method, a 2/3-rule for de-aliasing, 
and 512 3 collocation points on a [0, 2tt] 3 domain. The first DNS is of decaying turbulence which reaches a 
Taylor-microscale Reynolds number Re\ ~ 134 at the main peak of the enstrophy Z associated with the 
formation of the inertial sub-range. The Taylor-microscale A is defined in the usual way in terms of the 
energy spectrum E(k). The initial Fourier components of the velocity Uo(k) for the wave- vector k = \k\ 
are generated by applying random phases to the energy spectrum Eo(k) = Eok 4 exp{— 2k 2 }. The second 
DNS is a study of statistically steady turbulence which attains Re\ ~ 182 ; the forcing term f u (x,t) is 
specified most simply in terms of f u (k,t) whose spatial Fourier components are: 

f u (k,t) = v ^ f ~ k K (k, *), E u {k f , t)=j2 E u(k, *), (ii) 

^2E u (k f , t) k < kf 

where is the Heaviside function and kf = 2 is the wave number above which Fourier modes are not 
forced. This forcing term maintains a constant Grashof number Gr = L 3 V/v 2 = 4.9 x 10 T : for a similar 
forcing term that holds the energy injection fixed see Sahoo, Perlekar and Pandit (2011). 

For the decaying DNS, a small inertial subrange forms at t = 100 when the enstrophy Z reaches its 
main peak. Assuming E u (k) = Ko(k)ek~ 5 ^ 3 , the pre-factor Ko(k) is roughly 1.5 for about half a decade 
of wavenumbers. Similar to Fig. [l] (second panel), Fig. [2] (first panel) shows the time-dependence of 
the kinetic energy E, enstrophy Z and its skewness — S u . The second and third panels in Fig. [2] show 
Dm versus time t for m = 2, ... ,9 and = {w^ 1 ||k;||oo) °° with = \ for both the decaying and 
forced DNS calculations respectively. The second and third panels also show that D m < D m+ i and thus 
demonstrate the generality of Fig. [l]of ^TTJ 
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FIG. 3: Parameters of statistically stationary 
forced simulations : included are the resolution N, 
.Re a, the resolution parameter kmaxV an d the forc- 
ing type (see text). 
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FIG. 4: Scaling of the Q, m as a function of -Re^ for 
forced stationary isotropic turbulence with resolutions up 
to 4096 3 . Lines are for m — 1 (circles), 2 (squares), 3 (tri- 
angles), 4 (stars), 5 (left triangles), 6 (diamonds). Open 
and closed symbols correspond to EP and FEK forcing re- 
spectively. Dashed line is ~ -Re a ( see text). Note that 
for Re\ « 650 at 4096 3 with FEK forcing, moments up 
to fourth order (instead of sixth) are available from our 
database. 



IV. THE FOURTH SET OF SIMULATIONS : FORCED STATIONARY ISOTROPIC 

TURBULENCE 

The DNS data in this fourth set of simulations were obtained using a massively parallel pseudo-spectral 
code which achieves excellent performance on O(10 5 ) processors. The basic numerical scheme is that of 
Rogallo (1981). The time stepping is second-order Runge-Kutta and the viscous term is exactly treated 
via an integrating factor. Aliasing errors are carefully controlled by a combination of truncation and phase 
shifting techniques. The database includes simulations with resolutions up to 4096 3 and Taylor- Reynolds 
number up to Re\ « 1000 (Donzis et al 2012, Yeung et al. 2012). In order to maintain a stationary 
state, turbulence is forced numerically at the large scales. Since our objective is to assess the generality of 
the ordering of the moments D mi here, we show results using the stochastic forcing of Eswaran and Pope 
(1988) - denoted as EP - as well as a deterministic scheme described in Donzis and Yeung (2010) - denoted 
as FEK. In essence, this keeps the energy in the lowest wavenumbers fixed. For these two forcing schemes, 
the wavenumbers affected by forcing are confined to within a sphere k < hp, where hp is of order 2 or 3. In 
order to capture intense events, which are the main contributors to high-order moments, resolution issues 
have to be properly addressed. Motivated by the theoretical work of Yakhot and Sreenivasan (2004), 
resolution effects have been studied in Donzis (2012), and Yeung et al. (2012) with the conclusion that 
although high-order moments may be under-predicted using the standard resolution criterion - typically 
in simulations aimed at pushing up the Reynolds number - ratios of high-order moments are weakly 
affected by resolution issues. Small-scale resolution for a spectral simulation is typically quantified with 
the paramater k max r] where k max = ^/2N/3 is the highest resolvable wavenumber in a domain of size 
(27r) 3 with N 3 grid points. While the standard resolution is k max r] takes values between 1 and 2, results 
are presented from k max rj « 1.5 to 11, when available, which allows us to assess the effect of insufficient 
resolution. The Table in Fig. [3] summarizes those parameters of the DNS databased that have been used. 



A. The Dm— moments in forced stationary isotropic turbulence 

Even moments of vorticity Q m are shown in Fig. [4j As assured by Holder's inequality it can be seen 
that fL m +i > at all Reynolds numbers. The figure also shows the line ~ Re\ (dashed), which is 
the result of the dissipative anomaly. This is easily obtained from the kinematic relation (e) = v£t\ 
associated with isotropic turbulence and the well-known scaling (e) ~ Uq/L. It can then be shown that 
Qi ~ (Uo/^Re 1 / 2 ~ (Uo/L)Re\, where the well-known result Re\ ~ Re has been used. The DNS data 
in Fig. [4] agree with this scaling. As mentioned above, some resolution effects can be expected especially 
for high orders. Where data at nominally the same Reynolds number but different resolution is available, 
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FIG. 5: Scaling of moments of D m and ratios as a function of Re\ for forced stationary isotropic turbulence with resolutions up 
to 4096 3 . First panel : D m for m = 1 to 6. Second panel : Zoom of first panel to highlight the ordering of D m for m = 2 to 6. 
In both parts the dashed lines correspond to Re 6 x . Third panel : Ratio of moments D m _|_i/D m for m = 1 (squares), 2 (triangles), 
3 (stars), and 4 (left triangles) as a function of Re\. 



moments tend to be higher for higher values of k max r] (Donzis et al. 2008). This is clearer at higher 
Reynolds number (Re\ ~ 650 where two resolutions are available). Ratios of moments, however, are only 
weakly affected by resolution, which is also consistent with more recent results (Donzis et al. 2012, Yeung 
et al. 2012). In Fig. [5] the moments D m are shown as a function of Re\. For m = 1, one can again resort 

to using the dissipative anomaly with the definition D\ = (w^Vti) . The result is 

D\ = (L 2 i/{e)/^ 3 / 2 ) 2 ~ Re 3 ~ Re\ (12) 

which is seen in Fig. [5] (first panel). To see further details of higher order moments the second panel in 
Fig. [5] does not include D\. As in jpT| and j ]III| the data clearly shows the ordering D m +i < D m . The 
insensitive nature of moments to the type of forcing and the much weaker effect of resolution compared 
with Q m in Fig.[4]is also noted. The data also suggest that the ratio between successive moments decreases 
with m, which is consistent with the asymptotic behaviour of equation (J8|. This is more clearly seen in 
Fig. [5] (third panel) where the ratio of successive moments D m +i/D m is plotted for different values of m. 
Consistent with an ordering D m +i < D m , the ratio is always less than unity. As m increases, however, 
this ratio becomes increasingly closer to unity in agreement with Eq. [8] It is also interesting that these 
ratios appear to be independent of Reynolds numbers which suggest a regime I ordering with clustering 
of moments at high m also in the high-it^ limit. Resolution effects, while weak, can still be seen upon 
careful examination of the data, especially at high orders. However, for a given simulation, the ordering 
of regime I is unchanged with resolution. 



V. CONCLUDING REMARKS : THE DEPLETION OF NONLINEARITY 



The recent introduction of the Z) m -vorticity-moment-scaling (Gibbon 2011, 2012a, b), motivated by the 
time average Q, has suggested that they should be calculated through different numerical simulations. 
All four dats sets unexpectedly show that the Dm obey the ordering of regime I, namely D m +i < D m . 
This leads to the squeezing effect of ^ taking place such that ^ m+ i/^ m \ 1 and D m+ i/D m /l asm 



increases, which has an effect on the shapes of the PDF-tails, as remarked in <T C The ordering in the D m 
is strict although for m > 3, 4 the plots almost touch and replicate each other in shape as in Figs, [l] [2] 
and [5] even during intense events. It might be asked whether this is a viscous effect, or a strictly nonlinear 
effect, or the result of some surprising symbiosis between the two? Using a variation of the anti-parallel 
initial condition used in ^ TTJ new Euler calculations have repeated this observed ordering (Kerr 2012c), 
which implicates the nonlinear terms as the source. However, there is no evidence from Navier-Stokes 
analysis that such an ordering should hold, although no results exist that suggest it cannot. It is, of 
course, possible that a cross-over could occur between regimes I and II at Reynolds numbers higher than 
have been achieved in this work. 

Significantly D\ sits well above the other D m and does not appear to converge with them during the 
most intense periods : in Figs, [l] [2] and [5] D m lies on a log-scale with D 1 omitted. We are therefore 
justified in writing 

]nD m <a m ]nD 1 => D m < D\ m . (13) 
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Plots of a m for the first and second pair of simulations are shown in Fig. [6] Assuming a solution exists, 
the D m have been shown to obey (see Gibbon 2012a) 



D m < Dt 



D 



ra+1 



|m(4m+l) 



' ^2 < 



(14) 



n7oa m C2,m are algebraically increasing with m. 



where the c n ^ m within wi^ m = WQa m c 1 m and 

^2,? 

By dropping the negative term on the right hand side of (14), and replacing the D^-term with D m D\ G 
justified by (13), a time integration produces 



D m (i) < c m exp / D\ am dr < c m exp <( t 



j.l — 2arr 



t N 2a 

D x dr 



o 



2a m < 1 • 



(15) 



Fig. [6]shows that while there is a weak dependence of a m on both m and t, it nevertheless satisfies 2a m < 1 
in all cases. Because Leray's energy inequality insists that J* D\ dr < oo it is clear that the right hand side 
of ( [l5| is finite : any finite D m is sufficient for Navier-Stokes regularity. This regularization can thus be 
traced to the depletion of nonlinearity in (13) in regime I. Although regime II has not been observed, (14) 
shows that it is associated with time-decay of the D m . Specifically, if D m +i/D m > [ci )m C2, m ] 3 ^ 2m( ^ 4m ^^^ 
then D m < where [ci )m C2, m ] 3//2m ^ 4m+1 ^ \ 1 for large m. 
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